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Abstract 

Structural optimization of non-conservative systems with respect to stability 
criteria is a research area with important applications in fluid-structure in- 
teractions, friction-induced instabilities, and civil engineering. In contrast to 
optimization of conservative systems where rigorously proven optimal solu- 
tions in buckling problems have been found, for non-conservative optimiza- 
tion problems only numerically optimized designs were reported. The proof 
of optimality in the non-conservative optimization problems is a mathemat- 
ical challenge related to multiple eigenvalues, singularities on the stability 
domain, and non-convexity of the merit functional. We present a study of 
the optimal mass distribution in a classical Ziegler's pendulum where local 
and global extrema can be found explicitly. In particular, for the undamped 
case, the two maxima of the critical flutter load correspond to a vanishing 
mass either in a joint or at the free end of the pendulum; in the minimum, 
the ratio of the masses is equal to the ratio of the stiffness coefficients. The 
role of the singularities on the stability boundary in the optimization is high- 
lighted and extension to the damped well as to the case of higher 



Preprint submitted to 



January 4, 2011 



degrees of freedom is discussed. 
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1. Introduction 



Structural optimization of conservative and non-conservative systems with 
respect to stability criteria is a rapidly growing research area with important 



applications in industry 



2ll3| 



Optimization of conservative elastic systems such as a problem of optimal 
shape of a column against buckling is already non-trivial because some opti- 
mal solutions could be multi-modal and thus correspond to a multiple semi- 
simple eigenvalue which creates a conical singularity of the merit functional 
[sl. Despite these complications, a number of rigorous optimal solutions are 
known in conservative structural optimization. Nevertheless an increase in 
the critical divergence load given by the optimal design in such problems is 
usually not very large in comparison with the initial design 2|, |3|. 

In contrast to conservative systems, the non- conservative ones can loose 
stability both by divergence and by flutter. It is known that mass and stiff- 
ness modification can increase the critical flutter load by hundreds percent, 
which is an order of magnitude higher than typical gains achieved in opti- 

4|- 18|. For example, Ringertz jo] reported 



mization of conservative systems 



191] from 



an 838% increase of the critical flutter load for the Beck's column 
20.05 for a uniform design to 188.1 for an optimized shape. Recently Temis 
and Fedorov |l^ found for a free-free beam moving under the follower thrust 
an optimal design with the critical flutter load that exceeds that for a uniform 



2 



beam by 823%. We note that despite the very notion of the follower forces 



was debated in 
analogues 



20 



23 



24 



21 



22|, the Beck's column [l9[ as well as its discrete 



25| including the Ziegler's pendulum 26| remain popular 
models for investigation of mode-coupling instabilities in non-conservative 
systems and related optimization problems. 

In both conservative and non- conservative problems of structural opti- 
mization of slender structures, their optimal or optimized shapes often pos- 
sess places with small or even vanishing cross-sections. The known opti- 
mized shapes of the Beck's column or of a free-free rod moving under fol- 
lower thrust have almost vanishing cross-section, e.g., at the free end, which 
means vanishing mass of a finite element in the corresponding discretization 



0,0, 



11 



12 



17 
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Another intrigue of optimization of non-conservative systems is the 'wan- 
dering' critical frequency at the optimal critical load. During the optimiza- 
tion the eigenvalue branches experience numerous mutual overlappings and 



M 



27 



28|, 



29| with the tendency for the critical frequency to 



17 



This 



veerings 

increase and to correspond to higher modes [9|, [ij, [l5|, [ig, 
puzzling behavior of the critical frequency still waits for its explanation. 

In some problems, such as the optimal placement of the point mass along 
a uniform free-free rod moving under the follower thrust 0, 0], the local 
maxima were found to correspond to singularities of the flutter boundary such 
as cuspidal points where the multiple eigenvalues with the Jordan block exist 
12l |. In order for the last phenomenon to happen on needs at least three 



11 



modes [11|, [12] meaning that the two-modes approximations 



m 



are unable 



to detect such optima. This reflects a general question on model reduction 



3 



and validity of low- dimensional approximations in non-conservative problems 



discussed already by Bolotin [30] and Gasparini et al. [23| and recently 
raised up again in the context of friction-induced vibrations by Butlin and 
Woodhouse jsi l- 

The above mentioned phenomena make rigorous proofs of the optimality 
in the non-conservative optimization problems substantially more difficult 
than in conservative ones. To the best of our knowledge, the rigorously 
proven optimal solutions in optimization problems for distributed circulatory 
systems were not found. Although in the finite-dimensional case the situation 
is not much better, it looks reasonable to try to understand the nature of the 
observed difficulties of optimization on the final dimensional non-conservative 
systems that depend on a finite number of control parameters. 

Let us consider a circulatory system 

Mi + Kx = 0, (1) 

where dot indicates time differentiation, M is a real symmetric m x m mass 
matrix and K is a real non-symmetric mxm matrix of positional forces that 
include both potential and non-potential (circulatory) ones. The equation 
(II]) typically originates after linearization and discretization of stability prob- 
lems for structures under follower loads, in the problems of friction-induced 
vibrations and even in rotor dynamics when the damping is not taken into 
account 



30 



28 



32|. 



The characteristic equation for the circulatory system ([T]) is given by the 
modified Leverrier's algorithm jssl. In case when m = 2, it reads 

det MA^ + (trMtrK - tr(MK))A2 + det K = 0, (2) 



where A is an eigenvalue that determines stabihty of the trivial solution. 

Usually, the coefficients of the matrix K contain the loading parameter, 
say p, that one needs to increase by varying the coefficients of, e.g., the 
mass matrix. During this optimization process some masses can come close 
to zero so that the mass matrix can degenerate yielding det M = 0. As 
a consequence, some eigenvalues A can substantially increase. On the other 
hand such a singular perturbation of the characteristic polynomial may cause 
large values of the gradient of the critical load with respect to the mass 
distribution. We see that a problem of optimal mass distribution in a finite- 
dimensional circulatory system ([1]) in order to increase the critical flutter 
load, looks promising for explanation of the peculiarities of optimization of 
distributed non-conservative structures. However, it makes sense to tackle 
first not the most general system ([1]) but rather a particular finite-dimensional 
non-conservative system with m degrees of freedom. 



24 



n t 
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lis article, we propose to take an m-link Ziegler's pendulum 



23|, 



26| as a toy model for investigation of optimal mass and stiffness 



distributions that give an extremum to the critical flutter load. It appears 
that even the classical two-link Ziegler's pendulum had rarely been studied 
from this point of view in contrast to its continuous analogue — the Beck's 
column. The only example of such a study known to the author is contained 
in the book by Gajewski and Zyczkowski 

The paper is organized as follows. In the next Section we first consider op- 
timization of an undamped two-link Ziegler's pendulum. We find an explicit 
expression for the critical fiutter load as a function of the mass or stiffness 
distribution and demonstrate that in the space of the two mass coefficients 
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and the flutter load as well as in the space of the two stiffness coefficients 
and the ffutter load, the critical ffutter load forms a surface with a self- 
intersection and with the Whitney umbrella singular point. We consider the 
problem of optimal mass redistribution and find that the only two maxima of 
the critical fiutter load correspond to a vanishing mass either in a joint or at 
the free end of the pendulum; in the only minimum, the ratio of the masses 
is equal to the ratio of the stiffness coefficients. The maxima are attained at 
the singular cuspidal points of the stability domain and are characterized by 
a dramatic increase in the critical frequency of vibrations. Then, we write 
down the equations of motion of an m-link undamped Ziegler's pendulum 
and consider the case m = 3 in which we again find that the optimal mass 
distributions maximizing the critical flutter load correspond to vanishing of 
some of the three point masses. Other types of local extrema are also found 



that correspond to the points where the flutter boundary; 



las a vertical tan- 



11 



12l. or the points 



27 



28 



29j. Fi- 



gent such as cuspidal points with triple eigenvalues, cf. 
where the flutter boundary experiences a sharp turn, cf. 
nally we consider the problem of optimal mass distribution for a two-link 
Ziegler pendulum with dissipation. In conclusion, we formulate an optimiza- 
tion problem for an m-link Ziegler's pendulum and discuss some hypotheses 
on plausible optimal solutions and their expected properties. 



2. Structural optimization of the Ziegler's pendulum 

Let us consider the classical Ziegler's pendulum consisting of two light 
and rigid rods of equal length /. The pendulum is attached to a firm base- 
ment by a viscoelastic revolute joint with the stiffness coefficient ci and the 
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Figure 1: Undamped 2-link Ziegler's pendulum, (a) The critical load P(ci, C2) as a function 
of stiffness coefficients forms a conical surface in the (ci , C2 , P)-space (the case when mi = 
1712 — 1 and / = 1 is shown); (b) The critical load p{mi, 777,2) as a function of masses, forms 
a self-intersecting surface with the Whitney umbrella singularity at the point (0, 0, 2) of 
the (7771, ?7i2,p)-space (the case when ci = C2 = 1 is shown). 
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damping coefficient di. Another viscoelastic revolute joint with the stiffness 



32|. 



coefficient C2 and the damping coefficient d2 connects the two rods |26|, 
At the second revolute joint and at the free end of the second rod the point 
masses mi and m2 are located, respectively. The second rod is subjected to 



a tangential follower load P 



26 



32|. 



2.1. Undamped case 

Small deviations from the vertical equilibrium for the undamped Ziegler's 
pendulum are described by the equation (IB with the mass and stiffness 



matrices that have the following form 



26 



M = ;M '"M. K=r' + '^-^' I. (3) 

7712 m2 J y -C2 C2 

where x = {61,62)'^ is the vector consisting of small angle deviations from 
the vertical equilibrium position. 

Calculating the characteristic equation det(MA^+K) = for the Ziegler's 
pendulum without dissipation, we find 

mim2l^\'^ + {miC2 + 4m2C2 + Cim2 — 2Plm2)l^\^ + C\C2 = 0. (4) 
By direct calculation of the roots of the characteristic equation or by 



using the Gallina criterion 



24l | we find a critical surface that separates flutter 



instability and marginal stability domains 
(2/m2P - 4m2C2 - rriiCi - m2Ci)^ + (miCi - 171202^ = (miCi + 777-202)^. (5) 



The equation ([5]) defines a conical surface in the (ci, C2, P)-space when mi^2 
and / are fixed: flutter inside the cone. Fig. ^a). 
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However, in the (mi, m2, P)-space the critical surface (E]) looks differently 
and has a form of a self-intersecting surface with the Whitney umbrella sin- 
gularity at the (0, 0, 2c2//)-point. Indeed, expressing the critical load P from 
(jS]) we get 

^ _ 4m2C2 + (y^nlc^± ^/m^)^ ^ 2c2 , . 

~ 21^2 ~ ~r' 

Defining the new critical fiutter load as p := Pl/c2 we come to the more 
symmetric expression 

p = 2^UM±SX>2. (7) 



2 \ Y 777,2 V ^2, 

The case when Ci_^ C2 = 1, ttii = 2 and 7772 = 1 corresponds to the classical 



result of Ziegler 



26| 



p=l±V2. (8) 



The lower value of the critical load correspond to the boundary between 
marginal stability and fiutter while the higher critical load corresponds to 
the conventional transition from fiutter to divergence through the double 
zero eigenvalue with the Jordan block. 

The critical load ([7]) as a function of the masses p = p{mi, 7772) is plotted 
in Fig. mb). It is seen that the stability boundary has a self-intersection along 
a ray of the p-axis that starts at the a Whitney umbrella singularity with 
the coordinates (0, 0, 2) in the (777i,77i2,p)-space. Due to the symmetry of the 
expression ([7j), the critical load as a function of the stiffness coefficients, p = 
p(ci,C2), forms the identical surface in the (ci, C2,p)-space. In the following 
we will study the function p = p{mi, 7712). 
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Figure 2: Undamped 2-link Ziegier's pendulum with ci = 1, C2 — 1: The critical flutter 
load p{a) as a function of the azimuth angle a indicating the direction in the (rni,rn2)- 
plane. The point A is an absolute minimum of the flutter load: pa = 2, the point B 
corresponds to the local maximum: = 2 + ci/(2c2) with mi — 0, and the absolute 
maximum corresponds to a point C (not shown) with pc — +oo and m2 = 0. The point 
Z corresponds to the Ziegier's original design: mi/m2 = 2. 
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According to the inequality ([7]) the critical load is always not less than 
Po = 2. The minimum is reached when the masses satisfy the constraint 



miC2 = m2Ci. 



(9) 



Note that the equal stiffness coefficients Ci = C2 imply equal masses mi = m2. 
This situation corresponds to a uniformly distributed mass and stiffness in 



19|. 



continuous systems such as the Beck's column 

Usually, in the structural optimization problems the uniformly distributed 
stiffness and mass are considered as the initial design that is a starting point 
in optimization procedures. The critical load of the optimized structure is 
conventionally compared to that of the same structure with the uniform mass 



QflflflBlo, 
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13 



15 



17 



18|. 



and stiffness distributions 

Since p{mi^m2) is a ruled surface and thus p effectively depends on the 
mass ratio only, it is convenient to introduce the azimuth angle a by assuming 
rrii = cos a and m2 = sin a and to plot the critical load as a function of a. 
In Fig. [2] the curves p = p{a) bound the flutter domain shown in the light 
gray. When a tends to zero, which corresponds to the vanishing mass m2, the 
critical load increases to infinity. When a tends to ^ and, correspondingly, 
the mass mi is vanishing, then the critical flutter load increases to the value 



Pb 



Ici 
2 + — -. 

2C2 



(10) 



At the point B in the stability diagram of Fig. |2] the flutter boundary has 
a vertical tangent, which is a typical phenomenon in non-conservative opti- 



mization 



27 



28 



The lower part of the flutter boundary corresponds to the designs with a 
complex conjugate pair of pure imaginary double eigenvalues with the Jor- 
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dan block; the upper part is the designs with two real double eigenvalues of 
the same magnitude and different sign, each with the Jordan block. Above 
the upper flutter boundary lies the domain of statical instability or diver- 
gence with two unstable modes corresponding to two different positive real 
eigenvalues. At the point B the flutter boundary is tangent to the vertical 
part of the divergence boundary. To the right of this vertical line there is 
a pair of pure imaginary eigenvalues and a pair of real eigenvalues with the 
same magnitude and different signs. Transition from the stability boundary 
to the divergence boundary below the point B happens when a pair of pure 
imaginary eigenvalues goes out of the origin in the complex plane to infinity, 
merge there and return back along the real axis. This happens because at 
the boundary mi = 0, i.e. det M = 0. Similar divergence boundary exists at 
a = that corresponds to m2 = 0. Transition through the vertical line above 
the point B is accompanied by another eigenvalue bifurcation at infinity: two 
real eigenvalues of the same magnitude and different signs go out of the ori- 
gin in the complex plane in order to merge at infinity and then come back 
along the imaginary axis. The very point B corresponds to an antagonist of 
a quadruple zero eigenvalue with the Jordan block, i.e. to a quadruplet of 
complex eigenvalues that merge at infinity in the complex plane. 

To summarize, the popular initial design corresponding to uniformly dis- 
tributed mass and stiffness turns out to give an absolute minimum to the 
critical flutter load of the Ziegler's pendulum. The critical flutter load at- 
tains its local maximum, pb, for mi = at the singular cuspidal point B of 
the stability boundary where the flutter domain has a vertical tangent and 
touches the boundary of the divergence domain. Note that in \v}\ a local ex- 
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tremum of the flutter load for the free-free beam carrying a point mass was 
found to be at the cuspidal point on the flutter boundary too. The global 
maximum of the critical flutter load for the undamped Ziegler's pendulum is 
at infinity when 1712 = 0. 

The global maximum corresponds to a vanishing mass at the free end of 
the column which qualitatively is in agreement with the numerically founc 



optimized designs of the Beck's column available in the literature 



15 



16 



17 



13 



14 



isl . Indeed, all known optimized designs of the Beck's column 
are characterized by the vanishing cross-sections at the free end. Moreover, 
the gradients of the critical flutter load with respect to the mass or stiffness 
distribution of the Beck column are large, which is, again, in qualitative 
agreement with our stability diagram of Fig. [2J The most interesting is 
the fact that with the increase of the critical flutter load the higher and 
higher modes were reported to be involved into the coupling that indicates 



the onset of flutter j9|. 



13 



14 



15 



16 



17 



18| . Our simple model shows that 



this phenomenon seems to be natural for the optimal design that causes 
the degeneracy in the mass matrix that gives rise to the critical frequency 
that increases without bounds. It would be interesting to derive the optimal 
solutions for the optimization of the Ziegler's pendulum by means of the 
Pontryagin's maximum principle, which is used in numerical optimization of 



the distributed non-conservative systems, see e.g. 



10 



ll|. 



2.2. The m-link Ziegler's pendulum 

It would be very much desirable to extend our study to the case of the 
multiple-degrees-of-freedom Ziegler's pendulum. The corresponding models 
and recursive schemes for deriving the equations of motions where proposed 



13 



by Gasparini et al. [23j and by Lobas [25]. A three-link Ziegler's pendulum 



was considered by Gallina 



24|. 



We take the linearized equations of Lobas 
arbitrary masses and stiffnesses in the joints of an m- 
in contrast to the models of Gasparini and Gallina 



25| for its model deals with the 



m 



i Ziegler's pendulum 



23 



24j |. The mass and 



stiffness matrices of the Ziegler-Lobas model look like 



Em 
1=1 

Em 
i=2 



rrii 



E 
E 



rrii 



rrii 



Em 

i=m-l l^i 

Em V^*^ 

j=m-l l^i 



m \ 



rrii 



Em v^^*" - 

Em v-~^"^ 
^ i=m * /—ii=m * 



Em \— ^ 
i=m-l 1^ 



e: 



rrii 



E 



m 



rrii 
rrii 



(11) 



K 



Cl + C2 



PI -C2 
-C2 C2 + C3 - PI 













Cm— 1 ~l~ 






PI 

1 



Pl\ 
PI 

"Cm-1 + PI 
Cm 



For m = 2 the matrices (fTTj) and (fT^ are reduced to Note that 



det M = Z^"^ JJ 



rrii. 



[12) 



(13) 



i=l 



Let us consider the Ziegler-Lobas pendulum with m = 3 links. Now the 
mass at the free end is m^, while the masses m2 and mi are located at the 
third and second joints, respectively. The first joint connects the first rod 
with the basement. The length of each of the three rods is equal to /. The 
stiffness coefficients of the joints are C3, C2, and Ci, respectively. 
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Figure 3; Undamped 3-link Ziegler's pendulum with I — I, ci ~ C2 ^ ^ 1: The critical 
flutter load P{a) as a function of the azimuth angle a indicating the direction in the 
(rn2, ?7i3)-plane at the fixed radial distance r = 1 for (a) mi = 0, (b) mi — 10, and (c) 
TOi = 200 and for the fixed mi = 5 and (d) r = 1, (c) r = 0.65 and (f) r = 0.1. 
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For simplicity we assume that Ci = C2 = C3 = c. Then, the characteristic 
polynomial has the form 

aoA^ + aiA^ + a2A^ + 03 = 0, (14) 

with the coefficients 



rmim2?7i3, 
lit 



ao = 

01 = c/'*(6m2m3 + 5mim3 + mim2) — 2/^Pm3(mi + 7712), 

02 = WH'^m^i - 2(7m3 + m2)Pfc + (mi + 5m2 + 14m3)/^c^ 

03 = 



(15) 



Composing the discriminant matrix 

/ 



24| 








V 



(16) 



QjQ (X\ CI2 03 

3ao 2ai 02 
ao 0-1 02 03 
3ao 2ai 02 

Oq c^i ct2 c^s 
3ao 2ai 02 y 
and calculating the discriminant sequence consisting of the determinants of 
the three main minors of even order, we find that Ai = ?)l^'^m\m\m\ > 0. 
The expressions for determinants A2 and A3 = det S are rather involved 
and by this reason we omit them here. However, the numerical experiments 
evidence that the stability boundary is given by the equation A3 = for the 
stability condition A3 > implies A2 > 0. 

In Fig. |3] using the inequality A3 > we present the stability diagrams 
in the (a,P)-plane where the azimuth angle a in the (m2, m3)-plane is in- 
troduced by assuming m2 = rcosa and m^, = rsina. We assume the equal 



16 



lengths of the hnks / = 1 and the equal stiffness coefficients Ci = C2 = C3 = 1 
and vary the radial distance r in the (m2, m3)-plane and the mass mi. 

Since for m = 3 the critical surface P{m2, rris) is no more a ruled surface 
as it was in the case m = 2, the pictures in the {a, P)-plane change with the 
variation of the radial distance r that complicates the optimization problem. 
Nevertheless such diagrams are convenient for the analysis of the geometry 
of the stability boundary and thus for the identification of potential extrema. 
Moreover, the critical surface P{m2, m^) has a self-intersection along a ray of 
the P-axis that starts from the singularity Whitney umbrella as in the case 
of the 2-link pendulum. Therefore, at small values of m2 and rris the critical 
load can be locally approximated by a ruled surface. 

In the left column of Fig. [3] the radial distance r in the (m2, m3)-plane 
is fixed to r = 1 while the mass mi is increasing. As in the case m = 2, 
(marginal) stability is possible for a G [0,7r/2]. For mi = 0, two finite 
maximal values Pa and Pb are identified at a = {m^ = 0) and a = n /2 
(m2 = 0), respectively. Fig. M^a). Both maxima are attained at the cuspidal 
points of the stability boundary where it has the vertical tangents. However, 
the stability diagram changes when mi = 10, Fig.[3l^b). Again, local extrema 
exist at the boundary points a = {m^ = 0) and a = 7r/2 (m2 = 0), while 
the global maximum is at the point of the sharp turn of the flutter boundary 
with the vertical tangent near the cuspidal point Ci ~ (0.0403477, 11.961144) 
that corresponds to triple pure imaginary eigenvalues A ^ ±il. 1635243 with 



the Jordan block of third order, cf. 



11 



12| where at such a singular point 



a maximum of the critical flutter load was found for a free-free beam under 
the follower thrust. 
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With the further increase in the first mass up to mi = 200, the stabihty 
diagram converges to that similar to the diagram of the two-hnk pendulum, 
cf. Fig. [2] and Fig. |3t^c). This is not surprising because big inertia of the first 
joint makes the three-link pendulum effectively a two-link one. 

The right column in Fig. [3] corresponds to the fixed first mass mi = 5 and 
varied radial distance r in the (m2, m3)-plane. Small values of r correspond 
effectively to a two-link pendulum. That is why Fig. [3|^f) with r = 0.1 looks 
similar to Fig. ^c) and Fig. |2J 

Increase of r is accompanied by a complication of the stability diagram. 
In particular, two cusp point singularities originate corresponding to triple 
pure imaginary eigenvalues with the Jordan block. Fig. [3]J^d,e). Near these 
singularities, the stability boundary experience a sharp turn at the points Bi 
and B2 where the tangent to the boundary is vertical. At such points the 
eigencurves ImA(P) form crossing that can change either to avoided crossing 
or to the overlapping with the origination of a bubble of complex eigenvalues 
in dependence on the direction of variation of the azimuth angle a. This 



phenomenon was observed in many numerica 



optimization problems 



10 



13 



14 



analytically by Kirillov and Seyranian in 



11 



studies of non- conservative 



18| and had been described 



28| . Moreover, the critical load 



at these points can jump to a higher value corresponding to the merging 
of other (often higher) modes and is thus undefined ^, y, l7|. Nevertheless, 
these points could be local extrema of the merit functional, see 29| where 
the necessary conditions for that were derived. 

We stress that due to finite number of degrees of freedom and finite 
number of control parameters the stability boundary of an m-link Ziegler's 



18 



pendulum can be thoroughly analyzed both analytically and numerically. 
In particular, the coordinates of the singular points can be calculated with 
arbitrary precision and thus the issues of high sensitivity of the critical flutter 
load at the optima could be successfully resolved in this model, contrary to 
the optimization problems for distributed systems. A possibility to work 
with the singularities related to coalescence of more than two eigenvalues 
allows one to investigate qualitatively the question 'Should low-order models 



be believed?' 



3l| . Indeed, two-mode approximations work well far from such 



singularities, in their vicinity one has to take into account higher modes. 
2.3. Damped case 

In the presence of dissipation the equation of the Ziegler's pendulum is 



Mx + Dx + Kx = 



with the damping matrix 25|, |26| 



D 



di + d2 -d2 
-d2 d2 + ds 


















(17) 



■ ■ ■ djyi—i ~\~ djji d^_i 

\ ■ ■ • -dm-l dm J 

For the two-link damped Ziegler's pendulum with m = 2, we flnd the 
characteristic equation 



gqX^ + aiX^ + a2y? + OsA + 04 = 0. 



(19) 
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a) b) 




Figure 4: Damped 2-link Ziegler's pendulum with Z = 1, ci = C2 = 1: The critical flutter 
load P{a) as a function of the azimuth angle a indicating the direction in the (rni,rn2)- 
plane for (a) di = d2 = I and r = 1, (b) di = d2 = 1 and r = 0.4, (c) di = 1, d2 =0.1 
and r = 1, (d) di = 1, d2 = 0.1 and r = 0.1. 
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with the coefficients 



ai = P{mid2 + dim2 + 4m2(i2), 

a2 = (ii(i2 + mi/^C2 + 4m2/^C2 + Cim2/^ — 2P/^m2, 

03 = diC2 + Cid2, 

04 = C1C2. 



(20) 



Applying the Routh-Hurwitz criterion, we find the critical flutter load 

Aml{dlcl + dfcl) + did2{8m2{mi + 2m2)c2 + (cim2 — miC2Y) , 1 did2 



P 



For Ci = C2 



2(m2/(4m2(i2 + dim2 + mid2){cid2 + (iiC2) 



+ 



1, / = 1, mi = 2 and m2 = 1 it was found to be 
^"^ -r OOU1U2 -r -tu,2 _^ 

2' 



34| 



Adi + 33^1^2 + 4rf^ 1 , , 



2 7712/3 ■ 
(21) 



(22) 



2{<6d2 + di){d2 + di] 
The equation ( l22l) defines a surface with the Whitney umbrella singularity 
in the (rfi, 6/2, -P)-space which explains Ziegler's destabilization paradox by 



vanis hing dissipation 



1956 
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35| 



26| as it was first demonstrated by Bottema already in 
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M 



35 



we consider the 



In contrast to the previous studies, e.g. 
critical fiutter load ( pT]) as function of masses P = P{mi,m2) for the fixed 
damping distribution. 

In Fig. m the stability diagrams for the damped two-link Ziegler's pen- 
dulum are shown in the assumption of ci = C2 = 1, / = 1, mi = rcosa 
and 777-2 = rsina. The critical load P(mi, 777-2) does not constitute a ruled 
surface and thus the function P{a) depends on the radial distance r in the 
(ttii, 772-2)-plane. We see that the extrema again correspond to the boundary 
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points mi = and m2 = although the singularity at a = 7r/2 is an inter- 
section point. Nevertheless, with the increase of the number of degrees of 
freedom and control parameters new types of singularities will appear. The 
planar stability diagrams of Fig. H] depend on the damping distribution but 
do not tend to that of the undamped pendulum when damping goes to zero 



(destabilization paradox [32|, |35|, |36|). 



2.4- 'Problema novum ad cuius solutionem Mathematici invitantur' 

'A new problem that mathematicians are invited to solve' is a translated 
from Latin title of 1696 work by J. Bernoulli where he proposed the famous 
Brachistochrone problem js^. Supporting this good old tradition we would 
like to formulate the following optimization problem: 

Given a circulatory system (QP with the matrices M and K defined as 
in 177]) and ffM)- Find all local extrema, the absolute maximum of the 
critical flutter load P, and the corresponding extremal mass distributions 
{mi, 1712, ■ ■ -^rrim}- 

One can consider also the problem of optimal stiffness distribution or even 
vary both stiffness and mass. The same problems can be formulated for the 
damped pendulum with the damping matrix flTK]) . 

We expect that both in the undamped and in the damped case there exists 
a class of extrema, corresponding to the distributions {mi, m2, . . . , m^} with 
some masses mj = 0. It should be possible to find these optimal mass 
distributions explicitly perhaps with the use of the Pontryagin's maximum 
principle. It would be interesting to identify the singularities of the stability 
boundary that correspond to these extrema; some of them should be related 
to infinite eigenvalues A. 

22 



On the other hand, some local extrema should exist with the mass distri- 
butions that do not contain the vanishing masses rrij. It would be interesting 
to understand at which points — smooth or non-smooth — of the stability 
boundaries they are attained. Since the system is finite-dimensional and con- 
tains the finite number of control parameters with the clear physical meaning, 
the locations of the singularities corresponding to multiple eigenvalues can 
easily be found numerically with the high accuracy. In the vicinity of such 
points where at least three pure imaginary eigenvalues couple, the question 



3l| makes sense because here a one 



'Should low-order models be believed' 
more degree of freedom is crucial for the correct solution. 

Knowledge of the rigorously established optimal solutions at every m 
should shed light on the behavior of the optimal mass distributions and 
critical frequencies with the increase of the number of degrees of freedom. As 
a by-product, such a study will give an insight to the problem of dimension 



reduction and will serve a nice playground both for application o: 
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;he existing 
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40j and 



methods of nonsmooth analysis and eigenvalue optimization 
for their further development in the direction of more tight relation both 
with the singularity theory and the applications' needs. We expect that 
the proposed model optimization problem will yield useful recommendations 
for improvement of the algorithms of optimization of real non- conservative 
structures in industry. 

3. Conclusions 

In the present article we proposed to consider an m-link Ziegler-Lobas 
pendulum loaded by a follower force as a model problem for studying the 
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basic properties of optimization of circulatory systems with respect to sta- 
bility criteria. The model appears to be never studied from the point of view 
of structural optimization although its distributed analogue — the Beck's 
column — is a popular subject for investigation in numerical optimization of 
non-conservative systems. We have found the optimal solutions of the two- 
link undamped Ziegler's pendulum and established that they correspond to a 
vanishing mass either at the free end or at the second joint. At these designs 
the critical flutter load attains its maxima. A uniform design corresponds to 
an absolute minimum of the flutter load. These findings are qualitatively in 
agreement with the properties of known numerically optimized designs of the 
Beck's column. For the first time we have shown that the merit functional of 
the two-link undamped Ziegler's pendulum depending of two masses forms 
a singular surface with the Whitney umbrella singularity that governs opti- 
mization. We have formulated a problem of optimal mass distribution for 
an m-link Ziegler-Lobas pendulum and studied the case of m = 3 degrees 
of freedom. It turned out that in this case there exist extrema that corre- 
spond to vanishing masses at some of the joints. However, the new degree of 
freedom is responsible for the new types of extrema in the vicinity of the sin- 
gular points of the merit functional that correspond to triple pure imaginary 
eigenvalues with the Jordan block. Near such points, the stability boundary 
of the three-link pendulum qualitatively differs from that of a two-link one 
that clearly demonstrates the limits of a two-mode approximation. Finally, 
we have shown that although in the presence of damping the merit func- 
tional changes qualitatively, there still exist optimal solutions with some of 
the masses vanishing at the joints of the pendulum. 
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